Lecture 18

class: center, middle, inverse, title-slide

# Lecture 18
## Simple Linear Regression
### Psych 10 C
### University of California, Irvine
### 05/11/2022

---

## Simple linear regression

- We will start with another example using simple linear regression.

- We are interested in studying how the time that it takes participants to 
identify a geometrical object changes as a function of the degrees by which the 
object is rotated.

- For example, in a single trial of the experiment we present a participant with 
a geometric figure (left), and then after a a fixed period of time has passed 
(say 10 seconds) we present the participant another figure (right). The objective 
of the participant is to respond whether the figure is the same or not.

.pull-left[
.center[![](data:image/png;base64,#start.png)]
]

.pull-right[
.center[![](data:image/png;base64,#rotated.png)]
]

---

## Example: Mental rotation

- We presented participants with 6 figures rotated 10, 45, 90, 130, 
160 or 180 degrees respectively, and recorded the time it took them to identify 
if the second image was the same figure but rotated.

- We don't have access to a population of students to experiment with, 
so the age of participants in the experiment ranged from 7 to 21 years.

- The data from the experiment looks like this:

<div id="htmlwidget-001b18a9f5368894d9b2" style="width:100%;height:auto;" class="datatables html-widget"></div>
<script type="application/json" data-for="htmlwidget-001b18a9f5368894d9b2">{"x":{"filter":"none","vertical":false,"fillContainer":false,"data":[["1","2","3","4","5","6","7","8","9","10","11","12","13","14","15","16","17","18","19","20","21","22","23","24","25","26","27","28","29","30","31","32","33","34","35","36","37","38","39","40","41","42","43","44","45","46","47","48","49","50","51","52","53","54","55","56","57","58","59","60","61","62","63","64","65","66","67","68","69","70","71","72","73","74","75","76","77","78","79","80","81","82","83","84","85","86","87","88","89","90","91","92","93","94","95","96","97","98","99","100","101","102","103","104","105","106","107","108","109","110","111","112","113","114","115","116","117","118","119","120","121","122","123","124","125","126","127","128","129","130","131","132","133","134","135","136","137","138","139","140","141","142","143","144","145","146","147","148","149","150","151","152","153","154","155","156","157","158","159","160","161","162","163","164","165","166","167","168","169","170","171","172","173","174","175","176","177","178","179","180","181","182","183","184","185","186","187","188","189","190","191","192","193","194","195","196","197","198","199","200","201","202","203","204","205","206","207","208","209","210","211","212","213","214","215","216","217","218","219","220","221","222","223","224","225","226","227","228","229","230","231","232","233","234","235","236","237","238","239","240","241","242","243","244","245","246","247","248","249","250","251","252","253","254","255","256","257","258","259","260","261","262","263","264","265","266","267","268","269","270","271","272","273","274","275","276","277","278","279","280","281","282","283","284","285","286","287","288","289","290","291","292","293","294","295","296","297","298","299","300"],[473.09,990.5,1512.21,2188.01,2875.55,3351.05,345.17,919.23,1768.9,2278.36,3094.24,3140.63,647.95,980.59,1828.38,2361.32,2861.62,3169.47,823.36,1168.95,1985.51,2232.47,2964.1,3242.66,598.9,1282.95,1647.87,2346.88,3203.42,2865.62,568.26,1216.07,1816.12,2335.08,2903.23,3168.31,403.57,1231.51,2208.71,2074.05,2857.82,3375.42,380.97,835.12,1772.62,2408.47,2720.35,3008.75,674.41,985.46,1803.78,2565.8,3351.63,3661.35,535.13,1121.08,1845.48,2452.53,3213.2,2991.23,759.55,799.12,1413.94,2537.53,2584.57,3056.43,868.65,1114.68,1585.58,2410.4,2728.31,3143.79,275.81,964.79,2192.35,2376.24,2778.83,3096.57,557.65,887.14,1954.27,2807.35,2764.66,2835.23,640.6,1301.51,1497.88,2352.65,2956.36,3212.88,594.03,981.56,1470.4,2402.28,2918.86,3164.44,546.13,923.32,1911.51,2526.49,2905.81,3353.07,783.3,999.49,1862.38,2512.25,3070.52,3410.02,377.16,1186.13,1685.72,2486.5,2595.95,3267.89,538.22,959.13,1791.84,2649.41,2874.08,2913.84,646.7,921.52,2063.83,2026.91,2588.84,2821.4,505.7,1067.25,1788,2571.7,2821.3,3121.7,663.74,1282.44,1682.54,2327.84,2657,3204.03,911.57,1384.29,2161.29,2143.49,2868.69,2907.04,771.61,1086.15,1306.13,2592.39,3125.6,3206.01,1021.24,1185.16,1895.31,2411.05,2409.38,3585.79,830.02,718.4,1846.42,2424.88,3044.67,3341.43,709.36,1229.34,1795.07,2547.25,3101.87,3242.06,362.78,1182.59,1932.79,2716.41,2709.78,3187.17,428.94,1242.93,1844.94,2120.12,2853.21,3098.16,726.99,1112.06,1701.89,2448.48,3088.89,3109.54,1271.54,1494.67,1793.66,2706.98,2597.43,3035.37,503,1051.77,2019.42,2051.07,2472.2,3066.5,687.99,1387.24,1784.23,2673.71,3037.31,2991.79,282.52,1034.95,1904.45,2049.01,2668.67,3045.61,458.56,1200.22,1531.67,2391.09,2876.22,3120.63,612.07,1071.11,1472.93,2584.77,2962.07,3462.06,425.34,1038.85,1902.53,2345.63,3100.87,3483.46,690.09,1249.59,1824.97,2738.97,2397.99,3069.79,623.32,1367,1885.75,2167.33,2552.08,3287.96,515.84,1219.27,1959.88,2292.47,2824.6,3038.59,633.51,1610.78,1756.82,2320.61,2935.32,3345.57,730.82,1497.45,1671.93,2623.6,2939.25,2959.53,657.48,1475.4,1954.87,2530.61,3100.19,3489.75,771.53,804.1,2073.12,2618.09,2906.54,3394.56,316.12,1369.54,1955.4,2480.85,3016.8,3355.38,421.04,1159.88,1995.49,2369.21,3074.37,2980.35,397.24,643.73,1911.43,2478.08,3256.13,3150.92,1095.2,1341.51,1402.62,2294.35,2480.44,2677.38,492.69,825.98,2156.42,2476.23,3047.28,2966.63],[10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180,10,45,90,130,160,180],[19,19,19,19,19,19,18,18,18,18,18,18,20,20,20,20,20,20,13,13,13,13,13,13,21,21,21,21,21,21,15,15,15,15,15,15,11,11,11,11,11,11,17,17,17,17,17,17,13,13,13,13,13,13,9,9,9,9,9,9,11,11,11,11,11,11,7,7,7,7,7,7,14,14,14,14,14,14,12,12,12,12,12,12,21,21,21,21,21,21,17,17,17,17,17,17,12,12,12,12,12,12,10,10,10,10,10,10,12,12,12,12,12,12,11,11,11,11,11,11,19,19,19,19,19,19,15,15,15,15,15,15,19,19,19,19,19,19,15,15,15,15,15,15,17,17,17,17,17,17,16,16,16,16,16,16,13,13,13,13,13,13,9,9,9,9,9,9,15,15,15,15,15,15,19,19,19,19,19,19,10,10,10,10,10,10,8,8,8,8,8,8,17,17,17,17,17,17,7,7,7,7,7,7,19,19,19,19,19,19,14,14,14,14,14,14,10,10,10,10,10,10,8,8,8,8,8,8,17,17,17,17,17,17,7,7,7,7,7,7,15,15,15,15,15,15,7,7,7,7,7,7,16,16,16,16,16,16,8,8,8,8,8,8,11,11,11,11,11,11,16,16,16,16,16,16,8,8,8,8,8,8,20,20,20,20,20,20,20,20,20,20,20,20,15,15,15,15,15,15]],"container":"<table class=\"display\">\n  <thead>\n    <tr>\n      <th> <\/th>\n      <th>response_time<\/th>\n      <th>angle<\/th>\n      <th>age<\/th>\n    <\/tr>\n  <\/thead>\n<\/table>","options":{"pageLength":4,"rownames":false,"columnDefs":[{"className":"dt-right","targets":[1,2,3]},{"orderable":false,"targets":0}],"order":[],"autoWidth":false,"orderClasses":false,"lengthMenu":[4,10,25,50,100]}},"evals":[],"jsHooks":[]}</script>

---

## Example: Mental rotation

- Response times are saved in milliseconds (ms), thus the large numbers.

- For now we will assume that the age of the participant has no effect on 
response times.

- This is how we usually approach data analysis with linear models. Remember 
that our research question focuses on the effect of the angle of rotation on 
response time, not on the age of participants.

- Because we are interested in the effect of angle of rotation we will start 
with a simple linear regression.

---

## Null model

- The first model that we need is the Null model, which assumes that there is no
effect of rotation angle on the time it takes participants to identify that the 
figure has is the same but has been rotated (response time).

- Remember that the model is formalized as:
`$$y_i \sim \text{Normal}(\beta_0, \sigma_0^2)$$`

- Were `$y_i$` represents the *i-th* response time in the experiment, `$\beta_0$` is
the expected response time according to the Null model, and `$\sigma_0^2$` is the 
variance of response times according to the Null model.

- As we have mentioned before, the best guess for the expected value of a Normal 
distribution (when there is only one parameter) is the average of all response 
times in the experiment.

- Therefore, `$$\hat{\beta}_0 = \bar{y} = \frac{\sum_{i=1}^{n} y_i}{n}$$`

---

## Regression model

- The second model that we need is the simple linear regression. This model 
assumes that the expected response time changes as a function of the degrees
by which the figure has been rotated.

- In other words, the model assumes that the **angle of rotation** is a good 
predictor of response times.

- The model is formalized as follows:
`$$y_i \sim \text{Normal}(\beta_0+\beta_1x_i, \sigma_1^2)$$`

- Were `$y_i$` represents the *i-th* response time.

- `$\beta_0$` is the intercept or the response time to a figure that has been 
rotated by 0 degrees.

---

## Regression model

- `$\beta_1$` is the slope or the change in response time for a 1 degree increase 
the in the rotation of the figure, for example, the change in response time to 
a figure that has been rotated 10 degrees in comparison to one that has been 
rotated 11 degrees.

- `$x_i$` represents the angle of rotation of the *i-th* figure presented to the 
participant.

- Finally, `$\sigma_1^2$` represents the variability of response times with 
respect to their expectation `$\beta_0 + \beta_1x_i$`.

- Notice that this model assumes that the expected response time will be 
different depending on the angle of rotation of the figure, in other words, 
we expect response times to be around one value for figures rotated 10 degrees
and to have a different expectation when the figure is rotated 10.1 degrees.

- In other words, our prediction about the expectation of the response times 
are a continuous variable that follows a straight line!

---

## Regression model

- The prediction of the linear model for a given angle of rotation `$x_i$` are:
`$$\hat{\beta}_0 + \hat{\beta}_1x_i$$`

- We can obtain the estimators `$\hat{\beta}_0$` and `$\hat{\beta}_1$` using the 
`lm()` function in R.

- Remember that the function in R takes two arguments "**`formula =`**" which has 
to be written as: 
`$$\text{dependent-variable} \sim \text{independent-variable}$$`

- The second argument is "**`data =`**" which is the name of the object that 
contains our observations.

---

## Predictions: Null model

- Let's start by adding the predictions of the Null model to the data and 
calculating the squared error by observation.

```r
# Total sample size
n_total <- nrow(mental_rotation)

# Prediction of the Null model
null <- mental_rotation %>% 
  summarise("pred" = mean(response_time)) %>% 
  pull(pred)

# Add prediction to data and calculate error
mental_rotation <- mental_rotation %>% 
  mutate("prediction_null" = null,
         "error_null" = (response_time - prediction_null)^2)

# Calculate the Sum of Squared Error of the Null model
sse_null <- sum(mental_rotation$error_null)

# Calculate the Mean Squared Error of the Null model
mse_null <- 1/n_total * sse_null

# Calculate the BIC for the Null model
bic_null <- n_total * log(mse_null) + 1 * log(n_total)
```

---

## Predictions: Linear regression

- Now we can add the prediction and error of the "rotation" model.

```r
# Obtain estimators for beta0 and beta1
betas <- lm(formula = response_time ~ angle, data = mental_rotation)$coef

# Add prediction and error of the linear model
mental_rotation <- mental_rotation %>% 
  mutate("prediction_rotation" = betas[1] + betas[2] * angle, 
         "error_rotation" = (response_time - prediction_rotation)^2,
         "epsilon_rotation" = (response_time - prediction_rotation))

# Calculate the SSE of the linear model
sse_rotation <- sum(mental_rotation$error_rotation)

# Calculate the MSE of the linear model
mse_rotation <- 1/n_total * sse_rotation

# Calculate the variance accounted for by the linear model
r2_rotation <- (sse_null - sse_rotation)/sse_null

# Calculate the BIC of the linear model
bic_rotation <- n_total * log(mse_rotation) + 2 * log(n_total)
```

---

## Results

- Now that we have calculated all the values that we need we can summarize them
on a table in order to compare our two models and answer our research question.

| Model | Parameters | MSE | `$R^2$` | BIC |
|-------|:----------:|:---------------------:|:-----:|:--------------------:|
| Null  | 1          | `$8.8083\times 10^{5}$` | NA    | 4112|
| Angle of Rotation  | 2          | `$4.2866\times 10^{4}$` | 0.95 | 3211|

- The results indicate that the angle of rotation of a figure is a good predictor
of the time it takes to participants to identify if a rotated figure is the same 
as the one they saw at the start of the trial.

- We can also interpret the estimated values of the parameters in the model.

---

## Interpretation of model parameters

- When we interpret the parameters of a linear model we usually write something
like this:

- The estimated value of the intercept was 453, which means
  that the it takes participants 453 milliseconds on 
  **average** to identify a figure that **was not rotated** (rotated 0 degrees).
  
--

- The estimated value of the slope was 15, this suggests
  that it takes participants on **average** 15 milliseconds
  **more** to identify a figure **for each additional degree** that the figure 
  has been rotated by.
  
--

- Notice that when we interpret the parameters we always have to say 
"**on average**", this is because the model predictions are about the expected
response times, which means that there is error around those predicted values.

- Furthermore, the interpretation has to be made by making reference to the 
variables on the study, for example, "**expected response time in milliseconds**"
of "**for each additional degree of rotation**". This is because each estimate
is associated with a different variable in the study.

---

## Evaluating the linear model

- In the code of the linear model we added an additional calculation called 
**`epsilon_rotation`**, this variable is just the difference between observation 
and model prediction without squaring the value.

- We have said before that if we add the difference between observation and 
model prediction that difference will be equal to 0 by definition. This is
true regardless of the situation, however, we can use those differences to 
evaluate our model in a more specific way.

- In general we want those errors to have two properties:

1. The difference between observations and predictions follows approximately
  a Normal distribution centered at 0. This means that if we make a histogram
  of the difference, the histogram should be centered at 0 and be symmetric.
  
--

2. There are no patters, in the difference between observation and prediction 
  of the model. This means that regardless of how we define our *x-axis* on a 
  graph, the difference between observation and prediction should look like 
  white noise (there are no obvious patterns).
  
---

count: false
 
### Histogram of the difference between observation and prediction
.panel1-epsilon-hist-auto[

```r
* ggplot(data = mental_rotation)
```
]
 
.panel2-epsilon-hist-auto[
![](data:image/png;base64,#lec-19_files/figure-html/epsilon-hist_auto_01_output-1.png)
]

---
count: false
 
### Histogram of the difference between observation and prediction
.panel1-epsilon-hist-auto[

```r
  ggplot(data = mental_rotation) +
* aes(x = epsilon_rotation)
```
]
 
.panel2-epsilon-hist-auto[
![](data:image/png;base64,#lec-19_files/figure-html/epsilon-hist_auto_02_output-1.png)
]

---
count: false
 
### Histogram of the difference between observation and prediction
.panel1-epsilon-hist-auto[

```r
  ggplot(data = mental_rotation) +
  aes(x = epsilon_rotation) +
* geom_histogram(color = "#E72F52", fill = "#E72F52", binwidth = 40,
*                aes(y =..density..), )
```
]
 
.panel2-epsilon-hist-auto[
![](data:image/png;base64,#lec-19_files/figure-html/epsilon-hist_auto_03_output-1.png)
]

---
count: false
 
### Histogram of the difference between observation and prediction
.panel1-epsilon-hist-auto[

```r
  ggplot(data = mental_rotation) +
  aes(x = epsilon_rotation) +
  geom_histogram(color = "#E72F52", fill = "#E72F52", binwidth = 40,
                 aes(y =..density..), ) +
* theme_classic()
```
]
 
.panel2-epsilon-hist-auto[
![](data:image/png;base64,#lec-19_files/figure-html/epsilon-hist_auto_04_output-1.png)
]

---
count: false
 
### Histogram of the difference between observation and prediction
.panel1-epsilon-hist-auto[

```r
  ggplot(data = mental_rotation) +
  aes(x = epsilon_rotation) +
  geom_histogram(color = "#E72F52", fill = "#E72F52", binwidth = 40,
                 aes(y =..density..), ) +
  theme_classic() +
* xlab("Response time - Prediction")
```
]
 
.panel2-epsilon-hist-auto[
![](data:image/png;base64,#lec-19_files/figure-html/epsilon-hist_auto_05_output-1.png)
]

---
count: false
 
### Histogram of the difference between observation and prediction
.panel1-epsilon-hist-auto[

```r
  ggplot(data = mental_rotation) +
  aes(x = epsilon_rotation) +
  geom_histogram(color = "#E72F52", fill = "#E72F52", binwidth = 40,
                 aes(y =..density..), ) +
  theme_classic() +
  xlab("Response time - Prediction") +
* ylab("")
```
]
 
.panel2-epsilon-hist-auto[
![](data:image/png;base64,#lec-19_files/figure-html/epsilon-hist_auto_06_output-1.png)
]

---
count: false
 
### Histogram of the difference between observation and prediction
.panel1-epsilon-hist-auto[

```r
  ggplot(data = mental_rotation) +
  aes(x = epsilon_rotation) +
  geom_histogram(color = "#E72F52", fill = "#E72F52", binwidth = 40,
                 aes(y =..density..), ) +
  theme_classic() +
  xlab("Response time - Prediction") +
  ylab("") +
* theme(axis.title.x = element_text(size = 20),
*       axis.title.y = element_text(size = 20))
```
]
 
.panel2-epsilon-hist-auto[
![](data:image/png;base64,#lec-19_files/figure-html/epsilon-hist_auto_07_output-1.png)
]

---
count: false
 
### Histogram of the difference between observation and prediction
.panel1-epsilon-hist-auto[

```r
  ggplot(data = mental_rotation) +
  aes(x = epsilon_rotation) +
  geom_histogram(color = "#E72F52", fill = "#E72F52", binwidth = 40,
                 aes(y =..density..), ) +
  theme_classic() +
  xlab("Response time - Prediction") +
  ylab("") +
  theme(axis.title.x = element_text(size = 20),
        axis.title.y = element_text(size = 20)) +
* stat_function(fun = dnorm,
*               args = list(mean = 0, sd = sqrt(mse_rotation)))
```
]
 
.panel2-epsilon-hist-auto[
![](data:image/png;base64,#lec-19_files/figure-html/epsilon-hist_auto_08_output-1.png)
]

---

count: false
 
### Difference between observation and prediction
.panel1-epsilon-index-auto[

```r
* ggplot(data = mental_rotation)
```
]
 
.panel2-epsilon-index-auto[
![](data:image/png;base64,#lec-19_files/figure-html/epsilon-index_auto_01_output-1.png)
]